Spatiotemporal evolution in a (2+l)-dimensional chemotaxis model 
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Abstract 



Simulations are performed to investigate the nonlinear dynamics of a (2+l)-dimensional chemotaxis model of Keller-Segel (KS) 
type with a logistic growth term. Because of its ability to display auto-aggregation, the KS model has been widely used to simulate 
self-organization in many biological systems. We show that the corresponding dynamics may lead to a steady-state, divergence in 
^vq a finite time as well as the formation of spatiotemporal irregular patterns. The latter, in particular, appear to be chaotic in part of the 
^jPj range of bounded solutions, as demonstrated by the analysis of wavelet power spectra. Steady states are achieved with sufficiently 
^ large values of the chemotactic coefficient (x) and/or with growth rates r below a critical value re- For r > r^, the solutions of the 
differential equations of the model diverge in a finite time. We also report on the pattern formation regime for different values ofx, 
r and the diffusion coefficient D. 
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1. Introduction 

The Keller and Segel (KS) model fT) is one of the most ex- 
tensively studied chemotaxis models; it describes, for instance, 
the aggregation of cellular slime molds driven by chemical at- 
traction. Because it gives rise to "auto-aggregation", the KS 
model is used to describe many other phenomena including as- 
trophysical ones |l2l|3]|4l. Its mathematical properties, such as 
the existence of global solutions in time for certain choices of 
parameters and the finite time divergence for other choices, are 
understood to some extent (|5]|6l|2l. In particular, it is known 
that solutions of the one-dimensional case may diverge in a fi- 
nite time only if the chemoattractant does not diffiise, but in 
higher dimensional spaces divergences are common |7|. Also, 
Painter et al have shown that the KS model allows pattern for- 
mation and spatiotemporal chaos in one spatial dimension [6J, 
but these features have not been investigated in full detail in the 
(2H-l)-dimensional case. 

Therefore, in this paper, we consider the (2-i-l)-dimensional 
KS model with a logistic growth term, which in dimensionless 
form reads fS'jSI: 



du 
'dt 



V (DVm -xuVv) + ru(\- u). 



dv 
dt 



(1) 



(2) 
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where u{x, y, t) denotes the organism density and v{x, y, t) de- 
scribes the concentration of the chemical signal around the 
point (jc, y) at time t. The coefficients D and x respectively 
stand for the organism diffusion or motility and the chemotactic 
sensitivity. The validity of Eqs. ([T} and Q in the framework 
of chemotaxis is supported by some experiments on the Es- 
cherichia coli (E.coli) bacterium (see, e.g. Refs. Ii8,i9j) even if 
this model does not seem to reproduce all the observed chemo- 
tactic motions ^ . 

It is known that the solutions of Eqs. ([T]i and Q with cer- 
tain parameter choices reproduce an aggregation phenomenon 
which mimics chemotactic migration along the directions of 
gradients of self-produced chemicals. It remains to explore 
other regions of the parameters space in search of the param- 
eter values which lead to globally existing solutions, the finite 
time divergence as well as the formation of regular patterns. 
Furthermore, it is interesting to know whether chaotic patterns 
arise in regions of bounded solutions. 

In the present work, we investigate these questions by means 
of numerical simulations, and we show that the solutions of the 
KS equations in (2H-l)-dimensions, with a logistic growth term, 
diverge in a finite time when the cellular growth rate r exceeds 
a critical value r^. 

A detailed discussion of different versions of the KS model 
and its applications to chemotaxis can be found in the litera- 
ture, e.g., in Refs. B |6l |2l EO] El [H . Chemotaxis is a bio- 
logical process by which organisms alter their movements or 
orientations in response to chemical gradients and concerns, 
e.g., the behaviors of E.Coli, amoebae, endothelial and neu- 
trophil cells. Furthermore, it is of interest in phenomena such 
as cancer cell metastasis, embryogenesis, angiogenesis, tissue 
homeostasis, wound healing, immune response, progression of 
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Figure 1 : (Color online) The instability regions satisfying X { VD + V'') 
(positive values in the colorbar) are shown in the (D, r)-plane. 

diseases, finding food, forming the multicellular body of pro- 
tozoa etc. ifTSl fT4l [T5l . In brief, chemotaxis is the process by 
which many organisms probe their environment and determine 
their movements in response to what they have found. 

2. Nonlinear evolution 

Equations ([TJ and (|2]i have a nontrivial steady-state solution 
(m,v) = (1,1). A linear stability analysis |6| through which 
different parameter regimes can be found for the formation of 
spatiotemporal patterns shows that the corresponding modes are 

unstable for X > { + V'") with wave numbers satisfying 
ki < k < ki where 

- ^~2D~' ^ 2^ yl(D + r-xf-4Dr. (3) 

Thus, given a value of^, the stable (negative values) and the un- 
stable (positive values) regions can be represented in the (D, r) 
plane as in Fig. 1. In an interval [0, L] with homogeneous Neu- 
mann (zero-flux) boundary conditions we have k - nn/L with 
n = 0, 1, 2, whilst for periodic boundary conditions we have 
k = 2nn/L. We numerically investigate Eqs. ([T]l and (j2]) using 
zero-flux boundary conditions and with a time step dt = 0.001, 
grid size 200 x 200 in a simulation box of size x Ly where 
Lx = Ly = rm/k. We explore the evolution of the density u 
and of the chemical concentration v in a larger domain, with 
an initial condition (u,v) = (1,1 + /(e)), where /(e) < 0.01 is 
a spatially varying random perturbation. Furthermore, we per- 
form a wavelet analysis of spectra to show that chaotic patterns 
may form in some parameter region corresponding to bounded 
solutions. This is important for the excitation of many modes 
in the nonlinear dynamics. 

We note that for zero-flux boundary conditions on an inter- 
val of length Lx or Ly, the smallest non-vanishing mode is the 



Figure 2: (Color online) Numerical solutions of Eqs. ^ and |2| at ; = 200 
for two different values of r: (a) r = 1 and r = 1.5 with a fixed D = 0. 1 and 
X = (>■ The solutions first increase around the steady-state (1, 1) and then drop 
to a more or less constant value for a longer time. 

first, i.e. « = 1 and k = tt/L, = Jr/Ly. Then, using the condi- 
tion of instability, for both bounded and unbounded solutions, 
one may identify a parameter and a critical value below which 
patterns are regular, while they turn irregular if the parameter 
is higher than this critical value. If the solution is unbounded, 
it diverges in a finite time, as shown later For example, fixing 
the parameter values D = 0. 1 and = 6, we find that r has a 
critical value - 2, below which the wave pattern is regular or 
stable, and above which patterns exhibit irregular features. Fur- 
thermore, depending on the parameter values considered, there 
exist a critical value for ko - k/n, which dehmits the range 
of chaotic pattern formation. The latter also depends on the 
chemotactic nonlinearity for which the system exhibit finite so- 
lutions, as shown in the following subsection. 

To obtain steady-state oscillations, we consider the case with 
ko = 0.07 (i.e. n>l),D^ Q.\,x = 6, /ti = 0.45 < k < ka ^ 6.9 
and two different values of r: r - 1 (upper panel of Fig. 2) 
and r = 1.5 (lower panel of Fig. 2), which fall in the instability 
region (c.f. Fig. 1). Figure 2 shows that initially, the solution 
lies close to the steady-state value (1,1) before the organism 
kinetics dominates. Then, its density or the chemical concen- 
tration drops to a steady-state profile determined by the balance 
between the chemotactic nonlinear term and the dispersive dif- 
fusion term. We will later see that the bounded solution may no 
longer exist and finite-time divergence may develop if the value 
of r is increased above r - 2. 

2.L Variation of the domain length 

As the value of is lowered, i.e. as n increases in order that k 
falls in the interval {ki,k2), the domain length increases and the 
regular patterns disappear while irregular spatiotemporal struc- 
tures arise (see Fig. 3). Basically, the pattern selection leads to 
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Figure 4: (Color online) Normalized density wave i;,, , rescaled by the standard deviation of the sampling data (upper panel) and wavelet power spectra P,i (lower 
panel) as functions of the sampling time T in logarithmic scale corresponding to Fig. 4(d), indicating chaotic wave trains. 
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the excitation of many more modes involved in nonlinear inter- 
actions. As time advances, collision and fusion of patterns take 
place, giving rise to higher harmonic modes (merging) and to 
the emergence of new modes. We also note that increasing the 
domain size leads to higher concentrations of both u and v at 
a given time (t = 100 in Fig. 3), and enhances the disorder of 
patterns till a chaotic evolution is established for r < 2, the con- 
dition required for solutions to be bounded. This is illustrated 
by the analysis of wavelet power spectra portrayed as in Fig. 4. 

2.2. Variation of the growth rate 

As r exceeds its critical value rc = 2, the solutions quickly 
evolve and diverge in a finite time (Fig. 5), where the density 
and the chemical concentration are reported prior to the diver- 
gence time. Figure 5 shows four different states at different 
times, for D = 0.1, x - 6 and r - 4. We find that the time of 
divergence depends on the values of r. smaller r implies longer 
time before the solution diverges. For instance, the divergence 
time in Fig. 5 is f 34, whereas D = 0.1,;^^ = 6 and r - 3 im- 
ply a divergence at f s; 52. Further increase of the growth rate 
r leads to collision and fusion of patterns in a short time when 
the nonlinear logistic growth term dominates over the diffusion 
and chemotactic terms. One can calculate the exact divergence 
times for other sets of parameters. 

It is to be noted that since the system is two-dimensional in 
space, the energy transferred in the nonlinear interactions may 
not be due to the spatiotemporal chaos as in the previous sub- 
section, but could be primarily related to the wave collapse, i.e., 
to the divergence of solutions. However, this aspect needs fur- 
ther investigation, which will be subject of future work. 

2.3. Combinations ofD, x and r 

Next, we explore the pattern dynamics for a set of combi- 
nations of the system parameters, especially by increasing the 
diffusion as well as the chemotactic coefficient. Both param- 
eters control various distinctive features. Comparing the up- 
per panel of Fig. 6 with Fig. 3(d), we note that increasing x 
stabiUzes oscillations by reducing the values of u and v. On 
the other hand, increasing the diffiision coefficient, while keep- 
ing the other parameters fixed as in Fig. 3(d), one can observe 
chaotic patterns similar to those of Fig. 3(d) {c.f. lower panel of 
Fig. 6). This can be explained as follows. As increases, with 
D and r fixed, the region of instability determined by the values 
of k increases, since its lower limit k\ decreases while its upper 
limit k2 increases. In other words, n may be smaller or fewer 
modes may be involved in the nonlinear dynamics leading to 
chaos. On the other hand, when D increases with other param- 
eters fixed, the instability region is strongly reduced, since the 
opposite happens with the lower limit and upper limits for k, 
and n has to be larger in order for k to fall within the instabihty 
region k\ < k < ki. Thus, many more modes will participate in 
the collision and fusion process, and as a result spatiotemporal 
chaotic patterns may likely arise. Figure 7 represents the con- 
tour plot of the wavelet power spectra, corresponding to Fig. 
6, plotted against a normalized sampling time T (lower panel). 
This shows that the system is chaotic for a long time interval. 



It is also interesting to observe that the dynamics become more 
complex in nature at large values of t [c.f. Fig. 7(d)]. 

3. Conclusions 

We have explored the nonlinear dynamics of a (2-1-1)- 
dimensional chemotaxis model of KS-type. Specifically, by 
numerical simulations, we have focused on the existence of 
bounded solutions, solutions which diverge in finite -time as 
well as on the formation of different patterns including those 
arising from merging. 

In the parameter range, where the solutions remain perma- 
nently bounded, the existence of spatiotemporal chaotic states 
is demonstrated by the analysis of wavelet power spectra. These 
phenomena are important in the context of chemotactic migra- 
tion but so far, to the best of our knowledge, have not been ex- 
plored in two-dimensional KS models. We have shown that an 
increase of the chemotactic coefficient (< 10) and/ or keeping 
the growth rate (r < 2), stable solutions around the steady-state 
are generated after a long time {t > 200). Furthermore, tak- 
ing the growth rate beyond its critical value leads to the finite 
time divergence of the solutions for the density and the chemi- 
cal concentration. Conceming pattern formation and evolution 
in a large domain, suitable combinations of the parameters x, 
D and r, e.g. 5 < < 10 and 1 < r < 2, lead to bounded so- 
lutions which develop spatiotemporal chaos, whose complexity 
grows in time. 
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Figure 5: (Color online) Numerical solutions of Eqs. Q and J2j. Pattern formation and solutions diverging in a finite time are shown at different times before 
divergence: (a) t = 25, (b) t = 28, (c) ; = 30 and (d) t = 34. The parameter values are D = 0.\,x = 6 and r = 4. 
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Figure 6: (Color online) Pattern evolutions as solutions of Eqs. |TJ and J2j for different combination of parameters: (a) D = 0.1,;;' = 10 and /■ = 2 at ? = 100, (b) 
Z) = 1,^ = 6 and r = 2 at f = 100. 
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Figure 7; (Color online) Contour plot of the wavelet power spectra corresponding to Fig. 6 (lower panel) for (a) t=100, (b) t=500, (c)=1000, (d) t=5000. 
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